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Abstract 

The electronic structure of the zero-gap two-dimensional graphene has a charge neutrality point 
exactly at the Fermi level that limits the practical application of this material. There are several 
ways to modify the Fermi-level-region of graphene, e.g. adsorption of graphene on different sub- 
strates or different molecules on its surface. In all cases the so-called dispersion or van der Waals 
interactions can play a crucial role in the mechanism, which describes the modification of electronic 
structure of graphene. The adsorption of water on graphene is not very accurately reproduced in 
the standard density functional theory (DFT) calculations and highly-accurate quantum-chemical 
treatments are required. A possibility to apply wavefunction-based methods to extended systems 
is the use of local correlation schemes. The adsorption energies obtained in the present work by 
means of CCSD(T) are much higher in magnitude than the values calculated with standard DFT 
functional although they agree that physisorption is observed. The obtained results are compared 
with the values available in literature for binding of water on the graphene-like substrates. 
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Graphene, a two-dimensional layer of carbon atoms packed in a honeycomb lattice, is a 
unique physical object and it is under intent attention for the last several years due to its 
fascinating properties. pQ Starting from the first experimental works on the observation of 
ambipolar field effect and quantum Hall effect in graphene, [2] this material still continues to 
astonish scientists, demonstrating various interesting phenomena, like high carrier mobility, 
integer and half-integer quantum Hall effect, Klein tunneling, etc. |lj Furthermore, many 
interesting practical applications of this material have been proposed. Among them are high- 
frequency field-effect transistor, [3] flexible touch-screens, [1] single-molecule gas sensors [H] 
and many others. Almost all graphene-based devices, which are presently built or will be 
implemented in future, are based on the fact that in the neutral state of graphene the density 
of states at the Fermi level is zero and can be easily changed upon particular conditions, 
leading to the dramatic variation of the conductance response of the graphene conductive 
channel. Such conditions can be realized in different ways, like in the field-effect transistor 
via application of different-sign voltages to the back-gate electrode (controllable change of the 
conductivity of graphene from n- to p-type), adsorption of graphene on different substrates 
(e. g. graphene ribbon can connect two different metal contacts which induce different types 
of graphene doping, [6J thus allowing for formation oi n — p junctions in graphene), or 
via adsorption of atoms or molecules with different electron affinity on graphene (as was 
demonstrated in Ref. [5] graphene can be used for constructing precise gas sensors). In 
the last two cases the strength of the additional interactions on graphene, when adsorbed 
on a substrate or when adsorbing molecules, defines the changes in the density of states of 
graphene around the Fermi level and has to be carefully examined. 

Contributions to the binding of molecules on surfaces can be either chemical or physical 
in nature, i. e. chemisorption or physisorption, respectively. [7] Chemical binding typically 
implies a change in the electronic structure of both the molecule and the surface, either due 
to an ionic interaction, through charge transfer between substrate and adsorbate, or due to a 
covalent interaction, where orbitals deriving from the adsorbate and the substrate form new 
bonding and anti-bonding linear combinations. Physisorption, on the other hand, can arise 
through interaction of the permanent surface dipole with a permanent molecular dipole if it 
exists, through the interaction of the permanent dipole with an induced dipole or through 
interaction between fluctuating dipoles both in the adsorbate and the substrate. The latter 
contribution relates to the dispersion or van der Waals (vdW) interactions. Although, the 
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typical binding energy of physisorption is small (50 — 200 meV vs. more than 500 meV for 
chemisorption), this interaction plays an important role in nature and technique. |H] 

Contrary to the local chemical interaction, the dispersion interactions originate from 
long-range electron correlation effects and they are not captured by the standard density- 
functional theory (DFT) because of the local character of commonly used functionals. Con- 
sequently, DFT often fails to describe physisorption correctly. It is possible to improve 
the result when combining DFT with empirical forms for the van der Waals interaction |9] 
or when modifying existing functionals. [lOj A relatively high accuracy is accessible with 
a new exchange- correlation functional named van der Waals-density functional (vdW-DF), 
recently developed by Dion et al. [TT] However, taking into account the lack of systematic 
improvability within the DFT framework, a better way would be to employ methods beyond 
the DFT approach. 

For a correct and consistent treatment of physisorption interactions it is necessary to use 
high-level wave-function-based post-Hartree-Fock methods like the M0ller-Plesset perturba- 
tion theory [12j or the coupled-cluster method (CC). [13j One problem here is that a very 
accurate treatment, e.g. with the CC method, scales very unfavorably with the number of 
electrons in the system. From a physical point of view, however, this difficulty should be 
avoidable because the correlation hole around an electron is a fairly local object. For solids 
it is reasonable, therefore, to transform the extended Bloch orbitals of the periodic system to 
localized Wannier orbitals. The reformulation of the many-body wavefunction in localized 
orbitals defines the group of the so-called local correlation methods. One method of this 
type is the method of increments, originally proposed by StoU [H] and further developed by 
Paulus and coworkers (for Reviews, see Refs. [l5l[T6]). In this approach, a periodic HF calcu- 
lation is followed by a many-body expansion of the correlation energy, where the individual 
units of the expansion are either atoms or other domains of localized molecular orbitals. 
Another local approach to the correlations problem, firstly formulated for molecules, [T7j 
has been recently extended to periodic calculations. The latter has been implemented in 
a post-HF local-correlation computation code - CRYSCOR. Present time this code al- 
lows for inclusion of the correlation effects at the local M0ller-Plesset perturbation theory at 
second order (LMP2) level and only for non-conducting systems. Therefore, at the present 
level of progress, if coupled-cluster energies are desirable the quantum-chemical treatment 
for periodic systems is usually done using finite embedded clusters via application of the 
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method of increments, as is also done in the present work. At the same time, at the MP2 
level we essentially employ both techniques. 

Many theoretical studies have been focused on the investigation of physisorption on 
graphene-like substrates (for a Review, see Ref. jlS]), several of them are devoted to the 
H20/graphene system. [20H26] While ab initio data available in literature are limited to 
interaction energies calculated for polycyclic aromatic hydrocarbons and range from —104 
to — 249meV, DFT studies, employing periodic approach, often give controversial results 
(— 1.94eV [21] vs. from —18 to — 47meV |2Sl[2S])- This gives a hint, that the physisorption 
of water on graphene may not be reliably reproduced in DFT calculations and a more ac- 
curate quantum-chemical treatment is required. Here, we apply the CCSD(T) approach in 
the framework of the above method of increments to the adsorption of H2O on graphene. 

The structure of this paper is as follows: The next two sections describe two local corre- 
lation methods employed in our studies (Sec. [l]) and the computational procedure (Sec. 



Sec. Ill presents the main results. Conclusions are drawn in Sec. 



I. METHODS 



A. Method of increments for adsorption energies 

The method of increments combines Hartree-Fock (HF) calculations for periodic systems 
with correlation calculations for finite embedded clusters, and the total correlation energy 
per unit cell of a solid is written as a cumulant expansion in terms of contributions from 
localized orbital groups of increasing size. A detailed description of this approach can be 
found in Ref. [I5] In this section we outline briefly how this method can be applied for the 
calculation of adsorption energies. 

To quantify the molecule-surface interaction we define the adsorption energy as 

-E'ads = -E'cr+HaO " -E'er - E^^q = -E'fds + -^ad^) (1) 

where -EGr+H20 is the total energy of the Il20/graphene system, and Eq,. and Eh^q are the 
energies of the fragments at the same coordinates as in the H2O /graphene system corrected 
for the basis set superposition error (ESSE) according to the counterpoise scheme of Boys 
and Bernardi [27] . The Hartree-Fock energy, -E^g, is calculated for the periodic system in 
the standard way (for details, see Sec. [H]), and the -Ead" is calculated within the incremental 
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expansion. p8l 129] For latter quantity one has to take into account all orbital groups, that 
change due to the interaction of the molecule with the surface. For example, the correlation 
contribution to the interaction energy within the adsorbed molecule, the 1-body increment, 
can be defined as follows: rjA = e^'^'^^'^^ — e^^^^. As the total adsorption energy Ea^s in 
expression ([T|, also r]A and all r]i (the changes due to adsorption in the surface increments) 
terms are corrected for the BSSE. Other contributions occur in the system due to the 
simultaneous correlation of orbitals in groups from the molecule {A) and the surface 
e.g., 2-body increment rjAi = SAi — Sa — Si and 3-body increment rjAij = £Aij — ^a — — £j — 
VAi — VAj — Vij- -^adf "^^^ ^'^^ ^e Calculated as the sum of all r^-terms taken with the proper 
weight factors (according to their occurrence in the system under consideration): 

=VA + J2V^ + T.VA^ + T. V^J + E VMj + •.• • (2) 

i i i<j i<j 

All incremental calculations are performed with the program package MOLPRO [30] , 
using the coupled cluster treatment with single and double excitations and perturbative 
triples [CCSD(T)]. [HI] Further details are summarized in Sec. [IT} 

B. Local MP2 method 

In order to compare the results evaluated with the method of increments to the fully peri- 
odic model at the MP2 level, periodic LMP2 calculations have been performed additionally. 
The molecular local correlation scheme [1^ has been recently generalized to periodic systems 
and implemented in the CRYSCOR code. [18] It is based on the local representation of the 
occupied and virtual spaces by orthogonal localized orbitals (Wannier functions, WFs) and 
non-orthogonal projected atomic orbitals (PAOs), respectively. The WFs are constructed 
within the localization-symmetrization procedure and a posteriori symmetrized. [32] PAOs 
are generated in the reciprocal space by projecting the Fourier-images of atomic orbitals out 
of the occupied space. [33] 

Within the local approximation the virtual space for each orbital pair is restricted to 
so called orbital pair-domains, i.e. sets of PAOs, centered spatially close to either of the 
two WFs. Such a truncation of the virtual space is justified by the exponential decay (in 
case of insulators) of the doubles amplitudes with mutual WF-PAO separation. The list 
of the pairs is also truncated based on the decay of the correlation energy with inter- 
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orbital distance. In periodic systems the energy from the missing pairs can be extrapolated 
to infinity by fitting the corresponding Cg coefficients for orbital pairs of each type. [TH] 
The restriction of the virtual space according to the local scheme might lead to 1 — 2 % 
underestimation of the total correlation energy. Although this might amount to large value 
at the scale of relative energies, the local error is systematic and in most cases cancels almost 
completely, provided adequately large domains have been chosen. 

The zero band gap of graphene causes severe difficulties for the LMP2 method. The 
zero denominators lead to divergence of the perturbative MP2 estimate of the correlation 
energy in graphene. [51] Besides, the decay rate of the WFs is no longer exponential. How- 
ever, since in the present work we study the adsorption energy, the local scheme allows for 
partitioning of the latter into intra- graphene, intra-adsorbate and inter-graphene-adsorbate 
contributions. This makes possible to completely eliminate the problematic part, namely 
the divergent intra-graphene correlation, from the LMP2 treatment. The three contribu- 
tions to the interaction energy possess a clear physical interpretation. [SH] Inter-pair energy 
describes the dispersive effects, while the intra-contributions mainly show the reduction in 
the magnitude of the correlation energy due to the compression of the electronic densities 
caused by the exchange repulsion. For highly polarizable systems like graphene the inffuence 
of the latter in the correlation energy is expected to be very small, as is also conffrmed by 



incremental calculations (see Sec. III). Therefore the intra-graphene contribution can indeed 
be safely omitted. 

Finally, due to a poor localization of the tt-WFs of graphene the corresponding domains 
have to be taken sufficiently large. Test calculations showed that 30-atom domains for such 
WFs provide almost converged values for the inter-pair energies. For other specifications of 
the periodic LMP2 calculations see Sec. |TT] 

II. COMPUTATIONAL DETAILS 
A. Structural models 

To model a single H2O molecule on graphene, in periodic calculations (3 x 3) graphene 
supercells were used. When considering the correlation energy within the method of in- 
crements, a graphene sheet is mimicked by a finite fragment as shown in Fig. [T} All C-C 
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UP DOWN 




FIG. 1: The C58H2o-cluster chosen to model the graphene sheet when performing the correlation 
calculations via the method of increments (hydrogens are not shown). The H20-graphene arrange- 
ments are defined by the position of oxygen atom above the six-membered carbon ring (C/T/B) 
and the water orientation (UP/DOWN). 

distances are set to the experimental values, i. e. dc-c = o-Gv/ (2sm60°) = 1.421 A. The dan- 
gling bonds are saturated with hydrogen atoms and the C-H distances are set to 1.084 A. For 
the water geometry we used the following values: c/q-h = 0.9584 A and ^h-o-h = 104.45°. 
There was no geometry optimization performed in this work. For the water molecule, three 
adsorption sites are considered, namely, on top of a carbon atom (T), the center of a carbon 
hexagon (C), and the center of a carbon-carbon bond (B). For these positions, two different 
orientations of the molecule with respect to the graphene surface are examined, namely: the 
circumflex-like and the caron-like, denoted as UP and DOWN orientations, respectively (see 
Fig.g. 

B. Basis sets 

The basis sets employed are the polarized correlation-consistent valence-double-C basis 
(cc-pVDZ) [36] of Dunning for C and H of graphene and the aug-cc-pVTZ [36] basis for 
the water molecule. The only change needed to use these basis sets in the periodic case is 
the modification of p-type GTO from the original cc-pVDZ basis for the C atom: original 
Up = 0.1517 Bohr^^ was set to 0.17 Bohr~^. Apart from that, the standard molecular basis 
sets have been employed without any change. For test calculations we additionally employed 
a carbon triple- zeta-quality [6s3p2dlf] basis for graphene (to be denoted below as VTZ). 
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For this basis the s- and p-orbitals were taken from Ref. [37], with a modification that the 
2 smallest exponents for both s- and p- orbitals have been upscaled to 0.5 and 0.2 Bohr^^. 
The d- and /-orbitals in this basis are those of the cc-pVTZ basis set. |3S] 

C. HF and DFT calculations 

The periodic mean-field calculations for the studied system were performed with the pro- 
gram package CRYSTAL 09. [38j In order to obtain converged results for the HF binding 
energy, the default parameters were modified. For the pseudo-overlap tolerances IT0L4 and 
IT0L5, used for the prescreening of the exchange integrals, tighter values than usual have 
been employed: 15 and 80, respectively. The other CRYSTAL tolerances (ITOLl-3), used 
for the screening of Coulomb integrals, have been set to 7. The k mesh shrinking factors 
(isotropic Monkhorst net) have also been set to 12 corresponding to 74 A;-points to sample 
the irreducible BZ. The chosen values guarantee the stability of the HF and LMP2 solution 
for the considered system and basis sets. For compatibility, the DFT calculations reported 
here were performed with the same thresholds as used for the HF calculations. For the 
exchange-correlation energy functional we employed the generalized gradient approxima- 
tion as parameterized by Perdew et al. (PBE). [32] Some additional test calculations were 
performed by means of the VASP code. [40j 

In the context of DFT, the role of the long-range van der Waals interactions on the 
adsorption energy of the water-graphene system was considered by employing the semiem- 
pirical approach proposed by Grimme. [9j This method relies on corrections added to the 
DFT total energy and forces, based on a damped atom-pairwise potential CqR"^ (Cg rep- 
resents the dispersion coefficient for a given atom pair and R is the distance between the 
atoms) . 

In several cases, in order to obtain an improved value for the total energy, the vdW-DF 
functional [TT] has been applied to the charge density calculated by the VASP code. This 
was accomplished by utilizing the JuNoLo-code. [H] 
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D. Periodic LMP2 calculations and localization of the orbitals 

The periodic LMP2 calculations were done using the CRYSCOR code. [IH]. The orbital 
domains have been specified in the following way: for the WFs of the water molecules the 
whole molecular units were included in the domains, in graphene for the a- WFs the domains 
comprised 6 (first and second nearest neighbors), for the vr-WFs - of 30 atoms, respectively. 
Intra-graphene pairs were not included in the calculation. For the inter-graphene-water pairs 
the cutoff distance was set to 11 A. All the integrals have been evaluated via the density 
fitting approximation employing the direct space local fit and combined Poisson/Gaussian- 
type auxiliary basis sets of the quintuple-zeta quality. [12] The 20 x 20 fc-mesh was used for 
generating the PAOs. [33] 

When employing the method of increments, we use the localization procedure according 
Foster and Boys [42j as implemented in the program package MOLPRO. ^0] The localization 
of the orbitals in the graphene plane is somewhat delicate. The Foster-Boys procedure yields 
the a bonds between all atoms, and localizes the vr orbitals at each second C-C bond. This 
causes an interesting effect of making the rings formally inequivalent to each other from the 
point of the number of localized vr-orbitals belonging to them. Two possibilities can occur for 
the C-adsorption site, namely: it might correspond to (i) the center of the benzene-like ring 
(a combination of 3 cr and 3 tt localized orbitals (or "bonds" in the standard organic chemistry 
description); (ii) to the center of the ring, depleted of vr- localized orbitals (i.e. possessing a 
combination of rather 6cr-like localized orbitals). However, our test calculations show, that 
final contribution to the total correlation energy is not affected by this ambiguity, due to 
the relatively poor localization of the vr-orbitals and the weight-factor compensation. The 
same conclusion can be made when analyzing the results reported in Ref . [2^ , where for the 
interaction energy between H2S and graphene-like cluster of different sizes were used. 

III. RESULTS AND DISCUSSION 

In Fig. |2] we plot the dissociation curve for the C-UP geometry: the HF curve is purely 
repulsive, and the system is stabilized by electron correlation effects. With the method of 
increments one gets access to the individual contributions from the different orbital groups. 
In Fig. [3] the various contributions to the interaction energy of H20/graphene {d = 3.1 A) 
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TABLE I: Adsorption energies (i?ads) and equilibrium distances (do) computed for the wa- 
ter/graphene interface at different levels of theory, do is defined as a distance between graphene 
plane and either O or H, for UP and DOWN orientations, respectively. 



Structure 


CCSD(T) 




MP2 




PBE 


PBE-D2 




do (A) 


^ads (meV) 


do (k) 




do ik) 


EcAa fmeV) 


do (A) 


Eads (meV) 


C-UP 


3.06 


-108 


3.09 


-98 


3.69 


-20 


3.07 


-83 


C-DOWN 


2.61 


-123 


2.66 


-106 


3.52 


-19 


2.60 


-139 


B-UP 


3.05 


-102 


3.09 


-99 


3.70 


-18 


3.17 


-77 


B-DOWN 


2.64 


-118 


2.69 


-103 


3.68 


-18 


2.67 


-129 


T-UP 


3.06 


-110 


3.08 


-101 


3.70 


-19 


3.18 


-75 


T-DOWN 


2.69 


-135 


2.70 


-116 


3.67 


-19 


2.65 


-128 



are presented. For an estimate of the far-away contributions, which are neglected due to 
the cut-off the incremental expansion, we performed an C^i^'^-fit. Thus, the correlation- 
energy increments for the distances up to 7. 5 A were calculated explicitly, whereas the long- 
range contributions corresponding to distances up to = 12 A are obtained by the fitted 
extrapolation. The latter brings about 3 % to E^^J. While 2-body increments, implying the 
simultaneous correlation of orbitals in groups from the molecule and the surface {rjAi), yield 
the major part of the interaction energy, since they describe the vdW interaction between 
molecule and surface, the 1-body correlated contributions are found to be small (see also 
below). The latter observation is because changes in the localized orbitals are already 
captured to a large extent at the HF level. The 2-body increments of surface, rjij, and the 
investigated 3-body terms are even smaller, and can be neglected for further calculations. 
These data are in good agreement with the previously published result for the adsorption 
of H2S on graphene. [28] As in that case, the largest contributions to -Eads originate from 
the TT-orbitals closest to the adsorbed molecule, all a contributions account for only about 
one-quarter of the total adsorption energy. 

To prove that the chosen cluster is an adequate model to describe the studied system, the 
periodic LMP2 calculations were performed within the same basis sets as for the incremen- 
tal calculations. One observes good agreement between the periodic LMP2 and incremental 
MP2 correlation energies: the intra-pair interaction, that corresponds to J2VAi in terms of 
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FIG. 2: H20-graphene interaction energy as a function of substrate-adsorbate distance as obtained 
with different methods for the C-UP geometry (The periodic LMP2 calculation has been done for 
the distance, corresponding to the minimum of the potentail curve). 

the method of increments, obtained employing the periodic LMP2 amounts to — 162.3 meV 
while the incremental MP2 yields J2VAi = — 159.7meV (see Fig. [3]). A further contribution 
to compare is the intra-water interaction which is equal to 9.2 meV and 10.4 meV for the 
periodic LMP2 and incremental MP2 calculations, respectively. The intra-graphene contri- 
bution to the interaction is indeed insignificant (— 3.4meV). Therefore, the neglect of this 
term in the periodic calculations does not lead to a noticeable error. 

Our tests of basis set quality show that the use of the triple-zeta basis set for the graphene 
increases the attraction due to dispersion (— 162.3 meV vs. — 170.1 meV, for cc-pVDZ and 
VTZ basis sets, respectively) but at the same time increases repulsion of the HF contribution 
(57.7meV vs. 63.5 meV, respectively). Therefore, in total the increase in the interaction 
energy is insignificant. 

The resulting CCSD(T) as well as MP2 adsorption energies and equilibrium distances, 
as obtained for the six studied arrangements of the water molecule relative to the graphene 
layer, are listed in Tab. [Tj Evidently, the DOWN orientation is clearly more preferred in 
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this compared to the opposite one (i.e. UP) and the atop adsorption position is 

energetically most stable, although the variation in adsorption energy between different 
circumffex-like structures is not higher than 17meV. The general trend is the same when 
considering the MP2-results, albeit there is a deviation in i^ads up to 16 % as compared with 
the CCSD(T)-data. 

It is interesting to compare the presented CCSD(T) results with the ones obtained by 
means of DPT. As expected, the adsorption energies evaluated with the standard PBE func- 
tional are severely underestimated. The corresponding equilibrium distances are very large 
and energy the minima are shallow (see e.g. Pig. [2]). These results are in good agreement 
with previously published data. [25] The dispersion correction term represents the domi- 
nant contribution to the binding energy. When employing the PBE-D2 scheme one finds 
clearly observable energy minima at reduced equilibrium distances for all the considered 
geometries (Tab|l]). We note that the results obtained with GTO basis coincide with those 
evaluated with plane wave code. Whereas within the standard PBE approximation no en- 
ergetic preference regarding the adsorption site or orientation of the adsorbate has been 
observed, dispersion corrected DPT and the CCSD(T) results agree regarding the prefer- 
able orientation of the water molecule. Moreover, the calculated equilibrium distances are 
very similar to each other. Surprisingly, the energy difference between the circumfiex-like 
structure (DOWN-orientation) and its UP-counterpart is much lower when considering the 
CCSD(T) results. A further discrepancy between the data obtained by means of the method 
of increments and the PBE-D2 scheme is the preferable adsorption site: the position in the 
middle of the C-ring is shown to be the most stable one when using PBE-D2. To clarify 
the situation we performed single-point calculations for C-DOWN and T-DOWN geometries 
(PBE-D2 minima) employing the vdW-DP functional. As a result the T-DOWN geometry 
was found to be by 8meV more stable than C-DOWN, being in line with the CCSD(T) 
values. 

When considering previous estimates of water-graphene (graphite) interactions evaluated 
by means of highly accurate quantum-chemical methods, one notes, that they have been 
scattered in the wide range from —104 to — 249meV. The upper limit [21j is based on an ex- 
trapolation of the MP2 interaction energies calculated for polycyclic aromatic hydrocarbons 
with increasing size. This interaction energy is too large in comparison to the experimental 
results as well as numerous simulations (see e.g. [20] and references therein). It is known. 
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FIG. 3: The different contributions to the total interaction energy from the method of increments 
are plotted, for the C-UP arrangement of water molecule with respect to the graphene layer lattice 
when d = 3.1 A. CCSD(T) method is used when calculating the correlation-energy increments. 

that adsorption on graphene surface typically brings extra dispersion stabilization compared 
to polycyclic aromatic hydrocarbons of finite size because of an increased number of interact- 
ing carbon atoms. However, as has been shown in Ref. |2Dj, water is an exception from this 
rule: When going from coronene to graphene the interaction energy decreases, which can 
be explained by the fact that H2O is a molecule with a sizable permanent dipole moment. 
A combination of periodic DFT (PBE approximation) and CC calculations performed for 
water complexes with benzene, naphthalene and coronene yields water-graphene interaction 
energy equal to — 135meV and the equilibrium distance of about 2.68 A. Although these 
values seem to be in very good agreement with the result found in this work, the global 
minimum structure corresponds to the C-DOWN arrangement of the water molecule above 
the graphene surface. At the same time, the semi-empirical approach employed by McNa- 
mara et al. [43] predicts T-DOWN geometry to be the most stable one, when considering 
complexes of a single water molecule and a single-walled carbon nanotube, even though 
binding is slightly overestimated as compared to our result. 
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The experimental estimate for the adsorption energy of a single H2O molecule on graphite 
surface obtained by means of microcalorimetry is — 197meV. [H] Surprisingly, the interac- 
tion energies determined directly from the experimental contact angles of water droplets on 
the graphite surface are shown to range from —65 to — 97meV. |22] The results obtained in 
this work cannot be directly compared to these experimental values, first of all, since further 
graphite layers will affect the binding due to the long-range dispersion forces. Our test cal- 
culation for graphite modeled by two-layered structure show an increase of the adsorption 
energy by approximately 25 % (PBE-D2, C-DOWN structure). Furthermore, it may happen 
that the actual most stable arrangement of the H2O molecule above the graphite does not 
coincide with any of the six structures considered here. 



IV. CONCLUSIONS 



We performed CCSD(T) calculations by means of the method of increments for the 
adsorption of H2O on graphene. It has been shown, that the circumfiex-like orientation 
of water is more favourable than the caron-like one. Atop adsorption site is preferred by 
the water molecule and the most stable structure is characterized by an adsorption energy 
of — 135eV. Qualitatively this result is reproducible at MP2 level of theory, although the 
water-graphene interaction is systematically underestimated as compared to the benchmark. 
Both CCSD(T) and MP2 yield significantly larger adsorption energies than the previously 
reported periodic-DFT data. This is a consequence of the local nature of the commonly 
used functionals (e.g. PEE) and the result can be substantially improved when applying 
post-DFT dispersion corrected schemes. The semiempirical PBE-D2 treatment predicts 
reasonable adsorption energies and equilibrium distances, yet giving at the same time some 
discrepancies regarding the adsorption position compared to the CCSD(T) result. These 
discrepancies are eliminated when applying the vdW-DF functional on top of the charge 
densities calculated using the PBE approximation. From Tab. |T] it is clear that the water- 
graphene potential is very shallow, particularly so at the level of the CCSD(T) and MP2 
methods. Relative stabilities of the individual minima thus depend sensitively on the choice 
of method. Appreciating also the challenges of a proper water simulation |45j probably only 
a combination of accurate quantum-chemical calculations, DFT, and molecular dynamics 
can provide a reliable description of the water-graphene interface. 
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